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Abstract. These notes are based on lectures delivered by the authors at a Langeoog 
c/2 , seminar of SFB/TR12 Symmetries and universality in mesoscopic systems to a mixed 

audience of mathematicians and theoretical physicists. After a brief outline of the 
basic physical concepts of equilibrium and nonequilibrium states, the one-dimensional 
simple exclusion process is introduced as a paradigmatic nonequilibrium interacting 
particle system. The stationary measure on the ring is derived and the idea of the 
hydrodynamic limit is sketched. We then introduce the phenomenological Kardar- 
q I Parisi-Zhang (KPZ) equation and explain the associated universality conjecture for 

surface fluctuations in growth models. This is followed by a detailed exposition of 
\ a seminal paper of Johansson [59] that relates the current fluctuations of the totally 

' asymmetric simple exclusion process (TASEP) to the Tracy- Widom distribution of 

random matrix theory. The implications of this result are discussed within the 
framework of the KPZ conjecture. 

cn 
o. 

1. Introduction and outline 

o 

In statistical mechanics the study of systems that are far from equilibrium continues to 
attract considerable attention both in the physics and in the mathematics literature. As 
it turns out, exclusion processes first introduced by Spitzer |122] and their generalisations 
provide an excellent set of models that display rich and interesting nonequilibrium 
phenomena. Moreover, these processes are intimately related to a number of different 
models from statistical mechanics, combinatorics, probability and random matrix theory 
and a fruitful interplay between these fields, triggered by the seminal work of Baik, 
Deift and Johansson j6] , continues to produce spectacular results that yield very precise 
information on the behavior of such systems. 
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It is beyond the scope of this article to explain all or even a fair amount of these 
results in detail. Instead, we remain faithful to the original goal of the set of lectures on 
which this manuscript is based, namely to explain the fascinating developments in this 
field to an audience of scientists working in many different areas of mathematics and 
theoretical physics. We do not assume any significant acquaintance with the concepts 
of statistical mechanics or probability theory. In order to keep our presentation as 
elementary as possible we shall always focus on the simplest cases. Readers who wish 
to obtain more information on a specific topic will be referred to the literature. Here 
we will make use of the numerous reviews that have appeared recently in this area. 

In the remaining part of the introduction we provide a short outline of the topics 
that will be treated in this paper. 

In Sect. [2] we begin our discussion by explaining in general terms what physicists 
mean by the distinction between equilibrium and nonequilibrium systems, and by 
describing different types of nonequilibrium behavior. We then introduce in Sect. [3] 
a class of stochastic models, known as simple exclusion processes. They describe 
the stochastic motion of interacting particles on a lattice where the interaction 
is given by the exclusion property, i.e. two particles may not occupy the same 
site simultaneously. Interacting particle systems provide useful models for various 
nonequilibrium phenomena. Technically, they are (discrete or continuous time) Markov 
chains, and we will argue below in Sect. 13.21 that a simple yet precise criterion for the 
equilibrium vs. nonequilibrium character of a given system can be formulated within 
the general theory of Markov chains [144] . 

We then proceed to explain in more detail how one dimensional simple exclusion 
processes can be analyzed, specializing to the case where particles move only to 
neighboring lattice sites. If the probabilities to move right or left differ from each other 
one obtains the asymmetric simple exclusion process (ASEP). Mostly this paper will 
be concerned with ASEP and its subcase TASEP (totally asymmetric simple exclusion 
process) where all motion is unidirectional. The criterion of Sect. 13.21 identifies ASEP 
(and consequently TASEP) as nonequilibrium systems. In fact, ASEP has become a 
paradigmatic model for driven transport of a single conserved quantity and most of our 
discussion is focused on this class of models. For readers who are not familiar with ASEP 
and TASEP it might be useful at this point to have a look at the precise definition in 
Sect. O 

We begin our exposition of the analysis of ASEP with a discussion of stationary 
measures in Sect. 13. 3 1 where it is explained why the uniform distribution always provides 
a stationary measure in the simplest case of periodic boundary conditions. 

For a macroscopic description of the dynamics the notion of the hydrodynamic limit 
is introduced. The time evolution of the macroscopic particle density is then described 
by a hyperbolic PDE that can be solved by the method of characteristics (Sect. I3.4D . 
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After this has been established it is most natural to ask how much the process fluctuates 
around this macroscopic description. At this point it is useful to realize that ASEP is 
equivalent to a specific random model for surface growth which is known as 'corner 
growth' or 'single step' model (Sect. 13. 5p . Kardar, Parisi and Zhang (KPZ) conjectured 
in their seminal paper [63] that the fluctuation properties of a large class of (growth) 
models are universal. Assuming that KPZ universality can be applied to this particular 
growth model, one may obtain predictions for the scaling exponents of the fluctuations 
of the height of the surface and of the scaling exponent of the (spatial) correlation length. 
The KPZ conjecture and its implications for ASEP will be explained in Sect. HJ 

The KPZ conjecture is based on the KPZ equation which essentially adds stochastic 
driving to the hyperbolic PDE that describes the hydrodynamic limit. Unfortunately, 
the KPZ equation is difficult to analyze and it was only very recently that fluctuation 
results for solutions of the KPZ equation became available (cf. Sect. 18 .4f) . Before that 
an alternative approach to studying fluctuation properties was taken that turned out 
to be very fruitful. Rather than dealing with the KPZ equation itself one analyzes 
various specific models that are believed to belong to the KPZ universality class. A 
major breakthrough was achieved through a spectacular discovery on Ulam's problem 
for random permutations. 

In the early 60's Ulam raised the question to determine the asymptotics of the length 
of the longest increasing subsequence of a random permutation of the numbers 1, . . . , N, 
where it is assumed that all iV! permutations are equally likely. It took about 15 years 
to prove that the expected value of the length of the longest increasing subsequence 
behaves like 2y/N as N becomes large and there was strong numerical evidence that the 
fluctuations around the mean are of order N 1 ^ 6 . In a remarkable paper Baik, Deift, and 
Johansson [6] did not only prove that iV 1//6 was indeed the correct scaling, but they also 
identified the limiting distribution for the appropriately rescaled fluctuations. It came 
as a surprise that the limiting distribution coincides with the Tracy- Widom distribution 
that describes the fluctuations of the largest eigenvalue of matrices from the Gaussian 
Unitary Ensemble (GUE) as the matrix size tends to infinity. 

The article [6] became the starting signal for an explosion of research activities 
that continue until today. It became immediately clear that there is a variety of 
related combinatorial models (e.g. growth models, last passage percolation models, 
tilings, directed polymers in a random environment, tandem queues) that can be 
analyzed at the same level of detail and where Random Matrix distributions appear 
in the asymptotic description. For further information we refer the reader to the 
surveys pfl I129|. IBTR 168] 1125} 188] and the monograph [TT] . Two results that are of 
particular relevance for our discussion of KPZ universality were obtained independently 
by Prahofer and Spohn [99] and Johansson [59] . 

Prahofer and Spohn used the results in [6], and further developments in [TO], to 
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describe the height fluctuations of the polynuclear growth model (which is somewhat 
different from the corner growth/ single step model mentioned above, cf. Sect. E]). 
They obtained the scaling law predicted by KPZ universality. Moreover, they were 
able to identify the limiting distributions of the fluctuations. Again these are given by 
the Tracy- Widom distributions of Random Matrix Theory. A distinction needs to be 
made depending on the curvature of the surface. Flat surfaces lead to statistics from 
the Gaussian Orthogonal Ensemble (GOE), whereas the height-fluctuations of curved 
surfaces are described by GUE statistics. 

On the other hand, the results that Johansson presented in [59] immediately apply 
to the totally asymmetric simple exclusion process with step initial condition. One 
obtains the scaling exponent and the limiting distribution for the fluctuations of the 
particle flux, that are again given by the GUE Tracy- Widom distribution. We devote 
Sects. [5] -0 to explain [59J in great detail. This part of our presentation can be viewed 
as an expanded and more self-contained version of Sects. 3 and 4 of the review [110] by 
T. Sasamoto. 

According to the philosophy of our paper we explain the results of Johansson in the 
simplest case. More precisely, we consider the particle flux at the origin for a discrete 
time version of TASEP (dTASEP, cf. Sect. 13.11 (hi)) with step initial data. In Sect.[5]we 
formulate Johansson's result in Theorem 15 .41 and discuss its relation to KPZ universality. 

The proof of Theorem 15.41 naturally falls into two parts. The first part is of 
combinatorial nature. Via a representation by waiting times (Sect. 16. 1 p the problem 
is mapped to finding the longest subsequence in a list of alphabetically ordered two- 
letter random words that is weakly increasing in the second letter (Sect. 16. 2j) . By the 
Robinson-Schensted-Knuth algorithm (Sect. 16.31) one may represent the random words 
by pairs of Semi Standard Young Tableaux of the same shape (see Definition 16 . 2 j) . 
The advantage of this representation is twofold: On the one hand the length of the 
longest weakly increasing subsequence is simply given by the length of the first row of 
the corresponding Young Tableau. On the other hand there exist explicit formulae for 
counting the number of Semistandard Young Tableaux of a given shape, that can be 
derived using Schur polynomials (Sect. I6.4p . The result of all this reasoning is formula 
( l60l) . where A denotes the Vandermonde determinant (see also ( 1441) and Definitions 15.11 
andE2D. 

The second part of the proof of Theorem 15.41 is the asymptotic analysis of ( 160]) . 
The key observation is that the right hand side of (1601) has exactly the same structure 
as the formula for the distribution of the largest eigenvalue of GUE matrices. In 
particular the method of orthogonal polynomials (Sect. 17. II) can be applied to complete 
the proof of Theorem 15.41 in Sect. 17.21 The somewhat miraculous appearance of the 
Tracy- Widom distribution for the fluctuation of the particle flux of dTASEP is now 
explained on a technical level by the fact that Hermite polynomials (used for GUE) and 
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Meixner polynomials (used for dTASEP) look the same near their respective largest 
zeros after appropriate rescaling. The similarity of Hermite and Meixner polynomials is 
no coincidence. We briefly discuss the universal behavior of orthogonal polynomials in 
Sect. El 

As it was mentioned above the work of Johansson [59J and of Prahofer and Spohn 
[99] mark the beginning of a broad stream of research activities that continues to produce 
new and exciting results at a rapid pace. In Sect. [8] we briefly sketch and summarize 
those directions of recent research that are closely related to the question of KPZ- 
universality. The first generalization beyond [59J that we describe concerns the initial 
conditions. The results of Johansson apply for step initial conditions where every site 
to the left of the origin is occupied whereas every site to the right is empty. Based 
on their work on the polynuclear growth model Prahofer and Spohn |101] formulated a 
conjecture for the fluctuations of the flux for TASEP in the case of a general initial step 
profile with arbitrary constant particle densities pi and pr to the left resp. right of the 
origin. In Sect. I8.11 we explain this conjecture that has recently been fully established by 
Ben Arous and Corwin prj. Most remarkably, in a series of papers [130j- [T36] C. Tracy 
and H. Widom were able to extend some of these results to general ASEP. It should 
be pointed out that their proof is based on the Bethe Ansatz and does not use any of 
the nice but very special combinatorial identities that were crucial in the argument of 
Johansson. 

Furthermore we provide in Sect. [8] pointers to the recent literature regarding spatio- 
temporal correlations for (T)ASEP (Sect. 18.21) . interacting particle systems beyond 
(T)ASEP (Sect. 18. 3p . fluctuation results for the KPZ equation (Sect. 18. 4p . and physical 
experiments where KPZ behavior can be observed (Sect. 18.51) . We conclude the paper 
with a few remarks on integrability and universality. 

2. Equilibrium and nonequilibrium states 

The most fundamental concept of statistical physics is the distinction between 
microstates and macrostates in the description of systems with many degrees of freedom. 
To fix ideas, consider a classical iV-particle system (say, a gas in a box) described 
by a Hamilton function H(q,p) of position variables q = (<Zi, (Zdiv) an d momenta 
p = (pi, ...,pdN)- Particles move in a region Q C M. d of volume V = \Q\. Then a 
microstate is simply a point (q,p) in phase space, whereas a macrostate will be defined 
for the purposes of these lectures as a measure Px{q,p)dqdp parameterized by a set of 
macroscopic state variables (in short macrovariables) X. Here Px(l,p) is a function on 
phase space and dqdp denotes the canonical Liouville measure. 

Examples of macrovariables are energy, density, temperature or pressure. The 
macrovariables parametrizing the macrostate Px could have a dependence on space 
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and time, but to be useful they should be chosen such that they are slowly varying. 
This singles out in particular the conserved quantities of the underlying iV-particle 
system as candidates for macrovariables. The mapping from the microstate (q, p) to the 
macrovariables X is many-to-one, and the measure Px{q,p)dq dp gives the probability 
to find the system in a particular set of microstates (q,p) under the constraint that 
the macroscopic state is described by X. In principle, the time dependence (if any) of 
Px(q,p) is induced by the classical Hamiltonian dynamics of the microstate variables 
(q,p), but in practice well-chosen macrovariables are often found to satisfy autonomous 
evolution laws, such as the equations of hydrodynamics. The derivation of macroscopic 
evolution equations from microscopic Hamiltonian dynamics is the goal of kinetic theory. 
A (much simplified) version of this problem will be addressed below in Sect. 13.41 

In this perspective, equilibrium states are a subclass of macrostates which are 
attained at long times by a system that is isolated or in contact with a time-independent, 
spatially uniform environment. Characteristic properties of equilibrium states are that 

• the macrovariables X are time-independent and spatially homogeneous, and 

• there are no macroscopic currents (e.g., of mass or energy). 

The two most important examples of equilibrium states are the following: 

a. ) In an isolated system the energy E is conserved, the appropriate macrovariables are 

X = (E, V, N) and the equilibrium state is the measure induced by the Liouville 
measure on the energy shell {(q,p) '■ H(q,p) = E}. This is known in physics as the 
microcanonical measure. 

b. ) In a system at constant temperature T particles exchange energy with the walls 

of the box Q in such a way that the mean energy is fixed. The appropriate 
macrovariables are then X = (T, V, N) and the equilibrium state is of the form 

Pr tVtN ~ exp[-PH], P = l/T, 

known as the canonical measure. 

Having roughly characterized equilibrium states, we may say that nonequilibrium states 
arise whenever the conditions for the establishment of equilibrium are not fulfilled. As 
such, this definition is about as useful as it would be to define some area of biology as 
the study of non-elephants. We can be somewhat more precise by making a distinction 
between 

(i) Systems approaching equilibrium. By definition, the macrostate of such a system 
is time-dependent. In addition, systems in this class often become spatially 
inhomogeneous; an important and much studied case are systems undergoing phase 
separation [26J. 
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(ii) N on equilibrium stationary states (NESS). These systems are kept out of equilibrium 
by external influences. They are stationary, in the sense that macroscopic state 
variables are time-independent, and they may or may not be spatially homogeneous. 
In any case they are characterized by non-vanishing macroscopic currents. 

Examples for NESS are 

• Heat conduction. In a system with boundaries held at different temperatures there 
is a stationary energy current proportional to the temperature gradient (Fourier's 
law). 

• Diffusion. In a system coupled to particle reservoirs held at different densities there 
is a mass current proportional to the density gradient (Fick's law). 

• Electric conduction. Here particles are charged and move under the influence of 
a constant electric field. The particle current is proportional to the field strength 
(Ohm's law). 

Among these three examples, the first two can be further characterized as boundary 
driven, in the sense that the NESS is maintained by boundary conditions on the quantity 
that is being transported (heat, mass), whereas the last example illustrates a bulk-driven 
NESS maintained by an external field acting in the bulk of the system. 

NESS are the simplest examples of nonequilibrium states. Nevertheless, their 
description in the framework of classical Hamiltonian mechanics is conceptually subtle 
and technically demanding (see, e.g., |141j ). The main reason is that a Hamiltonian 
system under constant driving inevitably accumulates energy. In order to allow for the 
establishment of a steady state, dissipation has to be introduced through the coupling 
to an external reservoir, that is, a system with an infinite number of degrees of freedom. 

These difficulties can be avoided by starting from stochastic microscopic dynamics. 
While less realistic on the microscopic level, stochastic models provide a versatile 
framework for addressing fundamental questions associated with the behavior of many- 
particle systems far from equilibrium. The class of models of interest here are known 
in the probabilistic community as interacting particle systems. These are lattice models 
with a discrete (finite or infinite) set of states associated with each lattice site and local 
interactions. We focus specifically on exclusion processes, which are introduced in the 
next section. 

It is worth pointing out that the notion of equilibrium states in statistical physics, 
as outlined above, is much more restrictive than the usage of the corresponding term 
in most areas of mathematics, where an equilibrium is commonly understood to be 
any time-independent solution of some deterministic or stochastic time evolution. Thus 
NESS are equilibria in the mathematical sense. In Sect. 13.21 we will give a precise 
definition of what distinguishes physical equilibria from other time-independent states 
in the context of continuous time Markov chains. 
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3. An introduction to exclusion processes 

3.1. Definition 

The simple exclusion process was introduced in 1970 by Frank Spitzer |122j . Particles 
occupy the sites of a (^-dimensional lattice, which for the purposes of this discussion will 
be taken to be a finite subset Q C The particles are indistinguishable, which implies 
that a microstate or configuration of the system is given by 

V = Mxen e {0,l} n , 

where rj x = (1) if site x is vacant (occupied). The dynamics can be informally described 
as follows (for a detailed construction see |122[ 185]): 

• Each particle carries a clock which rings according to a Poisson process with unit 
rate (i.e., the waiting times between rings are exponentially distributed). 

• When the clock rings the particle at site x selects a target site y with probability 
Qxy{v) an d attempts to jump there. 

• The jump is performed if the target site is vacant and discarded otherwise; this 
step implements the exclusion interaction between particles and enforces the single 
occupancy constraint r) x — or 1. 

Together these rules define the exclusion process as a continuous time Markov chain on 
a finite state space; some general properties of such chains will be discussed in the next 
section. Interactions (beyond the exclusion interaction) can be introduced through the 
dependence of the jump matrix q xy on the configuration rj. Similarly, inhomogeneity 
associated with sites or particles can be introduced by letting the waiting times and the 
jump matrix depend explicitly on the particle positions or the particle labels, see [73] . 

We next restrict the discussion to the one-dimensional case with nearest neighbor 
hopping and without inhomogeneities or explicit interactions. Then 

qxy = q$y,x+l + (1 - q)$y,x-l- 

Informally, the particle attempts to jump to the right with probability q and to the left 
with probability 1 — q. The following cases are of interest: 

(i) q — 1/2 defines the symmetric simple exclusion process (SSEP). We will see below 
that this is really an equilibrium system. However, when defined on a finite lattice 
of sites x = 1, L and supplemented with boundary rates a, (3, 7, 5 which govern 
the injection (a, 5) and extraction (7, /3) of particles at the boundary sites i — 1 
and i = L, this model provides a nontrivial example for a boundary-driven NESS 
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(ii) q 7^ 1/2 defines the asymmetric simple exclusion process (ASEP). When considered 
on the one-dimensional ring (a lattice with periodic boundary conditions) the system 
attains a bulk-driven NESS in which there is a non-vanishing stationary mass 
current. This is the simplest realization of a driven diffusive system [116J. 

Note that the boundary conditions are crucial here. On a finite lattice with closed 
ends, which prevent particles from entering or leaving the system, an equilibrium 
state is established in which the bias in the jump probability is compensated by 
a density gradient; this is the discrete analog of a gas in a gravitational field, as 
described by the barometric formula. Another possibility is to consider a finite 
lattice with open ends at which particles are injected and extracted at specified 
rates [71] . This leads to a NESS with a surprisingly complex structure, see [20] for 
review. 

(iii) q = 1 (or 0) defines the totally asymmetric simple exclusion process (TASEP). In 
contrast to the case of general q, this process can also be formulated in discrete time 
[143] : In one time step t — > t + 1, all particles attempt to move to the right (say) 
simultaneously and independently with probability 7r e (0, 1]; moves to vacant sites 
are accepted and moves to occupied sites discarded. Such a discrete time dynamics 
cannot be defined for < q < 1, because it would lead to conflicts when different 
particles attempt to simultaneously access the same vacant site. 

For 7r — > the discrete time TASEP (dTASEP) reduces to the continuous time 
process in rescaled time nt, while for ir = 1 it becomes a deterministic cellular 
automaton which has number 184 in Wolfram's classification |142[ [77]. The case 
of general ir has been studied mostly in the context of vehicular traffic modeling 
[TT71I29]. 

Note that in terms of the waiting time picture sketched above, the discrete time 
dynamics corresponds to replacing the exponential waiting time distribution by a 
geometric distribution with support on integer times only. The exponential and 
geometric waiting time distributions are the only ones that encode Markovian 
dynamics [75]. The waiting time representation will play an important role in 
the exact solution of the dTASEP presented below in Sect]6l 



3.2. Continuous time Markov chains 

Before discussing some specific properties of exclusion processes, we outline the general 
setting of continuous time Markov chains (sec [106J for an introduction). Consider a 
Markov chain with a finite number of states i = 1, . . . , C and transition rates T^-. The 
rates define the dynamics in the following way: 

When the chain is in state i at time t, a transition to state j ^ i occurs in the time 
interval [t, t + dt] with probability T^dt. 
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The key quantity of interest is the transition probability 

Pki(t) = Prob [state % at t\ state k at 0] = Pi(t) 

where the initial state k is included through the initial condition -Pj(O) = 5^. The 
transition probability satisfies the evolution equation 

3& j¥=i 3 

which is known as the master equation in physics |138] and as the forward equation in 
the theory of stochastic processes [106J. Here the generator matrix 



Ay 



has been introduced. The master equation simply accounts for the balance of probability 
currents going in and out of each state of the Markov chain. To bring out this structure 
we rewrite ([I]) in the form 

= Kij ; Kij = ^jiPj ~ ^ijPi, (2) 
3 

where Ky is the net probability current between states i and j |144j . If the chain is 
irreducible, in the sense that every state can be reached from every other state through 
a connected path of nonzero transition rates, the solution of (CQ) approaches at long 
times a unique, stationary invariant measure P* determined by the condition 

J>^* = 0. (3) 

3 

The invariant measure is the left eigenvector of the generator matrix, with eigenvalue 
zero. Based on (T5j) we can rewrite (J3J) as 

J>* =0 with K*t — TjiP* — TijP*. (4) 

3 

Two classes of Markov chains may now be distinguished depending on how the 
stationarity condition (j3J) is realized: 

(i) K*j = V In this case the probability currents cancel between any two states 

W = ^iP*, (5) 

a condition that is known in physics as detailed balance. In the mathematical 
literature Markov chains with this property are called reversible, because (jSJ) implies 
that the weight of any trajectory (with respect to the invariant measure) is equal 
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to that of its image under time-reversal |1U6[ [55]. Detailed balance or, equivalently, 
reversibility is a fundamental property that any stochastic model of a physical 
system in equilibrium must satisfy, because equilibrium states are distinguished by 
invariance under time reversal 

(ii) K*a 7^ at least for some pairs of states i,j. Such a Markov chain is irreversible 
and describes a system in a NESS. 

Examples for both kinds of situations will be encountered in the next section. 
3. 3. Stationary measure of the exclusion process 

We consider the ASEP on a ring of L sites with a fixed number N of particles. The 
total number of microstates rj is then C = ( N ) and the transition rates are 

( q : (... •o...) -)■ (... o....) 

T(ri^rf) = l l-q : (... o •...) - > t ... • o...) (6) 
[ : else. 

Here (... • o...) denotes a local configuration with an occupied site (•) to the left of a 
vacant site (o), and it is understood that only configurations 77,7/ that differ by the 
exchange of a single particle-vacancy pair are connected through nonzero transition 
rates. The stationary measure P*{rf) is determined by the condition 

Y / T( V '^r ] )P*(r ] ') = J2 T (v^v')P*(v) Vt/. (7) 

rf rf 

As the simplest possibility, let us assume that the invariant measure is uniform on the 

state space, 

P*(r/)=(J^) =► K*( V , V ') = [T( V ' ^ V )-T( V ^ V ')] 

We discuss separately the symmetric and the asymmetric process. 

• q = 1/2 (SSEP). Here the rate q — 1 — q = 1/2 for all allowed processes, and 
for each allowed process the reverse process occurs at the same rate. We conclude 
that detailed balance holds in this case, K* = 0, and the SSEP is reversible as 
announced previously. 

• q 7^ 1/2 (ASEP). Because for any allowed process with rate q the reverse process 
occurs at rate 1 — q 7^ q and vice versa, detailed balance is manifestly broken, 
K* 7^ 0, and we are dealing with an irreversible NESS. However, we now show that 
the uniform measure (jHJ) is nevertheless invariant. To see this, consider the total 
transition rates for all processes leading into or out of a given configuration rj. We 
have 

= £r(i/ 77) = q KM + (i - qWM 

v' 
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where N *{rj) denotes the number of pairs of sites with a particle to the right of a 
vacancy in the configuration rj. Similarly 

itSffa) = J2 T ^ -+ if) = s^fa) + - ^-M- 

A little thought reveals that Af, (v) — A4»(^) for any configuration 77. Hence 
Ttotl 7 ?) = r°ot (77) for any g, and the stationarity condition (JTj) is satisfied for the 
uniform measure (jSJ). 

A few remarks are in order. 

(i) The invariance of the uniform measure (jHJ), and the fact that it is independent of 
the bias q, relies crucially on the ring geometry With open boundary conditions 
allowing for the injection and extraction of particles both the SSEP and the ASEP 
display nontrivial invariant measures characterized by long-ranged correlations and 
the possibility of phase transitions [38| [20] . For example, for the SSEP with 
boundary densities pi at x = 1 and pa at x = L one finds a linear mean density 
profile, as expected from Fick's law, but in addition there are long-ranged density- 
density correlations on the scale L, which take the form [38], 1123] 

HvitVL?) - Efo^EfeO = - ~^'\ pl - PR?- 

Here £, £' G [0, 1] are scaled position variables with £ < 

(ii) The invariant measure of the dTASEP on the ring is not uniform. Rather, one finds 
a Gibbs measure with repulsive nearest-neighbor interactions between the particles 
|143[ I115[ 1117] . This means that the probability of a configuration can be written 
as a product of pair probabilities, 

P* p (v)~l[p P (Vx,Vx+i), (9) 

X 

where the limit N, L — > 00 at fixed density 

p = E(jfc) = N/L (10) 

is implied and 

1 — \/l — 4np(l — p) , , 

Pp (o, 1) = Pp (i, o) = — v 27T py ^, (11) 

p p (0, 0) = 1 - p - p p (l, 0), p p (l, l) = p- p p (l, 0). (12) 

For 7r — > this reduces to a Bernoulli measure of independent particles (see 
Sect. I3.4p . whereas for 7r — >• 1 we have p p (l,0) — > (1 — |1 — 2p|)/2, which implies that 
p p (l, 1) for p < 1/2 and p p (0, 0) — > for p > 1/2. At n = 1 and mean density 
p = 1/2 the measure is concentrated on the two configurations r] x = [1 ± ( — l) x ]/2. 



Interacting particle systems, KPZ universality, and random matrices 



13 



(iii) The invariance of the uniform measure for the ASEP on the ring is an example of 
pairwise balance |119j . a property that generalizes the detailed balance condition 
(jSJ) into the form 

T(r] -> rj')P*(rj) = T(rj" -> rj)P*(rj"). 

This means that for each configuration rj' contributing to the outflux of probability 
out of the state rj there is a configuration rj" whose influx contribution precisely 
cancels that outflux. In other words, the terms in the sums on the two sides of ([71) 
cancel pairwise. 



S.Jf.. Hydrodynamics 

An important goal in the study of stochastic interacting particle systems is to understand 
how deterministic evolution equations emerge from the stochastic microscopic dynamics 
on large scales |124[ |66| [TBj. This is similar to the (much harder) problem of 
deriving hydrodynamic equations from the Newtonian dynamics of molecules in a gas 
or a fluid. The mathematical procedure involved in the derivation of macroscopic 
evolution equations for systems with conserved quantitities is therefore referred to as 
the hydrodynamic limit. Here we give a heuristic sketch of hydrodynamics for the ASEP. 

The key input going into the hydrodynamic theory is the relationship between the 
particle density p and the stationary particle current J. The particle current is defined 
as the net number of particles jumping from a site x to the neighboring site x + 1 per 
unit time, which is independent of x in the stationary state. From the definition of the 
ASEP we have 

J = qE[ Vx (l - Tfc+i)] - (1 - qMVx+i(l - rj x )] 

where expectations are taken with respect to the invariant measure. Since all 
configurations of iV particles on the lattice of L sites are equally probable, 

m N(L-N) 
- %+i)] = - Tfe)] = J L _ 1 ■ 

This is just the probability of finding a filled site next to a vacant site, which is obtained 
by first placing one out of N particles in one of L sites, and then placing one out of 
L — N vacancies in one of the remaining L — 1 sites. We conclude that 

where the particle density (TlOl) is kept fixed. Similarly 

= JJ^J ^P 2 = HVxMVy] for L -> oo (13) 
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for any pair of sites x 7^ y. This implies that in the invariant measure on the ring, 
restricted to a fixed finite number of sites, for L — > 00 each site is occupied independently 
with probability p (Bernoulli measure). 

We can now formulate the basic idea of the hydrodynamic limit [83]. Suppose 
that we start the ASEP at time t = from a Bernoulli measure with a slowly varying 
density p(x,0). Here "slowly varying" means that variations occur on a scale i ^> 1. 
Since the invariant measure of the ASEP is a Bernoulli measure of constant density, it 
is plausible that, if £ is chosen large enough, the evolving measure will remain close to 
a Bernoulli measure with a time and space dependent density p(x, t) at all times; and 
because the particle density is locally conserved, the evolution equation for p(x, t) must 
be of conservation type, 

J^ (M) + J^' (M) = °- (14) 
In the limit £ — > 00 we may expect, in the spirit of a law of large numbers, that the 
local particle current j(x, t) converges to the stationary current associated with the local 
density p(x, t), 

j(x,t) -»■ J(j>(x,t)), 

such that (114j) becomes an autonomous, deterministic hyperbolic conservation law 

t + ! J <^° < 15 > 

for the density profile p(x,t). Equation (ITS"]) is known as the Euler equation for the 
ASEP, because similarly to the Euler equation in fluid mechanics it lacks a second order 
"viscosity" term vd 2 pjdx 1 . It must be emphasized that such a term does not appear to 
leading order, when the hydrodynamic limit is carried out at fixed q ^ 1/2. It is only 
present in the weakly asymmetric case, which implies that q — > 1/2 in the limit £ — > 00 
such that £(q — 1/2) is kept fixed [37] . 

The Euler equation ( JT5|) has been rigorously established for a wide range of models, 
including cases in which the invariant measure and the current- density relation J(p) 
are not explicitly known [120] . We conclude this section by a brief discussion of the 
properties of the nonlinear PDE (ITS]) , assuming a general (but convex) current-density 
relation with J(0) = J(l) = 0. This includes in particular the dTASEP for which 

JO) = irp p (l,0) = i[l - v/l-47rp(l-p)]. (16) 



Shock formation. Hyperbolic conservation laws of the form ([15]) can generally be 
solved by the method of characteristics. To this end we first rewrite ( TTo']) in the 
form 

dp , ,dp n 
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where 



. . dJ 

c{p) = Tp 



A characteristic is a trajectory of a point of constant density, and the key 
observation is that the characteristics of ( ITT)) are straight lines. Denoting by v po (x, t) 
the position of a point of density po = p{x, 0) at time t, we have to satisfy the 
condition 

p(v Po (x,t),t) = po = p(x,0) 

at all times. Taking the time derivative of this relation and using (j!7p we see that 
the solution is 

v Po (x,t) = x + c(po)t, 

i.e. points of constant density travel at the kinematic wave speed (j!8j) . 
The convexity of the current- density relation implies that c(p) is a decreasing 
function of the density. As a consequence characteristics collide in regions of 
increasing initial density, dp(x,0)/dx > 0, leading to the formation of density 
discontinuities (shocks) in finite time. At this point the description by the PDE 
( fTBI) breaks down, but the speed V of a shock separating regions of density Pl on 
the left and pn > pl on the right is easily inferred from mass conservation to be 
given by 

v Apr)-J(pl) ^ (19) 

PR- PL 

Note that V — > c for pi — > p^. On the microscopic level shocks are represented 
by the shock measures of the ASEP [17J H5] . These are inhomogeneous invariant 
measures on Z which approach Bernoulli measures with density pi and pn for 
x — > — oo and x — > oo, respectively. The microscopic structure of shocks has been 
studied in considerable detail [39|. 

Rarefaction waves. If the initial density profile is a step function 

^-{: : <-» 

with pl > pr, a diverging fan of characteristics forms leading to a broadening, 
self-similar density profile 

x < c(p L )t 

x > c(p R )t (21) 
c(pl) < x/t < c(p R ), 

where the shape function </>(£) can be computed from the current-density relation 
J(p). Inserting the ansatz ( l21i) into ( fTTj) we see that 

m=c-\0- (22) 



p(x,t) 
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Figure 1. (Color online) Schematic of the mapping between a configuration of the 
ASEP and the corresponding height configuration h x+ i (bold blue line). Initially 
all sites x < in the ASEP are occupied, and all sites x > are vacant (step initial 
condition) . The ASEP occupation variables determine the height differences according 
to (|23|) . At the same time, the height increment h x+ i(t) — h x+ i(0) counts the net 
number of particles that have crossed the bond (x, x + 1) from left to right up to time 
t. In the figure two particles have crossed the origin, and hi = 2. Rotating the height 
configuration by 45° provides a representation of the net number of jumps a given 
particle has undergone as a function of the particle label, counted backwards from 
x = (the Vth particle is the one that was located at x = — I at time t = 0). 

For the continuous time ASEP the interpolating shape is linear, 




3.5. Mapping to a growth model 

The representation of the one- dimensional ASEP as a growth model seems to have been 
formulated first by Rost [107J. It is illustrated in Fig. [1] for a step initial condition, 
where all sites to the left of the origin are occupied (r) x = 1 for x < 0) and all sites to 
the right of the origin are empty (j] x = for x > 0). This initial condition will also play 
a central role below in SectJSJ 

The mapping assigns to every configuration r] = {i] x } of the ASEP a configuration 
of height variable {h x+ i}, where the shift of the index indicates that the height variable 
h x+ i lives on the bond connecting the ASEP sites x and x + 1. After fixing the height at 
a reference point, e.g. by setting hi = 0, the height configuration is uniquely determined 
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by the relation 

h x+ i -h x _i=--r) x . (23) 

The ASEP occupation variable encode the local slopes of the height profile, which take 
the values h x+ i — h x _i = ±|, hence the name "single step model" in the physics 
literature [921 M\. 

The step initial condition corresponds to the initial height profile 

h x+ i{0) = ^\x\, 

which accounts for the designation as a "corner growth model" in the mathematical 
literature. It can be seen from Fig. [1] that a particle jumping across a bond to the 
right (left) increases (decreases) the corresponding height variable by one unit. Thus 
the height h x+ i is an odd (even) multiple of \ for odd (even) values of x, and the height 
change h x+ i (t) — h x+ i (0) is equal to the net number of particles that have jumped across 
the bond (x, x + 1) from left to right up to time t (jumps from right to left are counted 
with a negative sign). Finally, for the special case of the step initial condition, the net 
number of jumps (forward minus backward) performed by a given particle can also be 
read off from the height profile (see Fig. [T]). 

The mapping ( l23l) is clearly not restricted to the step initial condition. Of particular 
interest are translationally invariant initial conditions, which can be constructed 
deterministically or stochastically. For example, to generate a deterministic initial 
condition of density p = 1/n, one simply places a particle at every n'th site of the 
lattice, and a stochastic initial condition is obtained by occupying sites independently 
with probability p. The two types of initial conditions differ in the roughness of 
the corresponding height configuration, which is quantified by the height difference 
correlation function 

G(r) = E[(Vr - h y f] - E[h y+r - h y ] 2 . (24) 

An ensemble of height configurations on Z is said to be smooth if lim^oo G(r) < oo and 
rough otherwise [79]. The deterministic initial conditions described above are smooth in 
this sense, whereas for the stochastic initial condition a simple computation using ( 1231) 
and (Tl3|) shows that 

G(r) = p(l-p)|r|. (25) 
4. The KPZ conjecture 

The asymmetric exclusion process and the equivalent growth model introduced in 
the preceding subsection are representatives of a large class of models, which was 
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brought to the forefront of research in nonequilibrium statistical physics in 1986 by 
a seminal paper of Kardar, Parisi and Zhang (KPZ) [63]. Working in the framework 
of a phenomenological stochastic continuum description, they formulated what may 
be called a universality hypothesis encompassing the fluctuation properties of a large 
class of different microscopic models^]. The classic period of research in this area has 
been extensively reviewed in the literature [791 ESI [72] . Here we aim to give a concise 
and simple presentation of the KPZ conjecture, in order to place the more recent 
developments (to be elaborated in the following sections) into their proper context. 

4-1. The Kardar- Parisi- Zhang equation 

We start from the hydrodynamic equation ( TT5|) with a general current-density relation 
J(p). Since we are interested in fluctuations around a state of constant mean density p, 
we write p(x, t) = p + u(x, t) and expand to second order in u, which yields 

du , .du , du , . 

a = " A "a? < 26 > 

where 

A-g(fl. (27) 

The linear drift term on the right hand side can be eliminated by a Galilei transformation 
x — > x — ct, which leaves us with what is known (for A = 1) as the inviscid Burger 
equation. 

Now fluctuations are introduced (in the spirit of fluctuating hydrodynamics [124]) 
by adding a random force to the right hand side of ( |26|) . In order to guarantee mass 
conservation, this term must take the form of a derivative —d(/dx of a stochastic process 
((x, t) in space and time. This is assumed to be a stationary Gaussian process with zero 
mean and a covariance function 

E[((x, t)((x', t')} = a- l a- l G[(x - x')/a x , (t - t')/a t ] (28) 

which vanishes beyond a small correlation length a x and a short correlation time a t . 
Usually one takes formalfyfjj] a X} a t — > 0, which reduces the right hand side of ( 1281) to a 
product of ^-functions, 

E[C(x, t)C(x', t')} -> D5(x - x')5(t - t') (29) 

and turns the process ((x, t) into spatio-temporal white noise of strength D. This rather 
violent driving has to be compensated by a small viscosity term ud 2 u/dx 2 with v > 0. 
Putting all ingredients together we thus arrive at the stochastic Burgers equation 

du d 2 u du d( d . . 

dt dx 2 U dx dx dx^ X ' ' 

| For an introduction to the idea of universality from a mathematical perspective see [34j . 

§ In the hydrodynamic context [51] it is also of interest to consider the solutions of (I5U)) on scales small 

compared to the spatial driving scale a x - 
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first introduced in the context of randomly stirred fluids [5T] and subsequently applied 
to fluctuations in the exclusion process by van Beijeren, Kutner and Spohn [137J. 

To establish the connection to growth models we proceed in analogy to the discrete 
case discussed in Sect. 13.51 We introduce the height function h(x,t) through the time- 
integrated particle current, 

h{x,t) — I j{x,s)ds, (31) 
Jo 

Supplementing this with the initial condition u(x, 0) = 0, it follows from the conservation 
law for u that 

dh 



* ™ + (33) 



and therefore 



dt dx 2 2 \dx 

which is precisely the KPZ-equation [63]. In general there is also a constant term on the 
right hand side of (133]) which has been set to zero. The defining feature of the equation 
is the quadratic nonlinearity on the right hand side, which is present whenever A 7^ 0, 
that is, when the current is a (generic) nonlinear function of the densitjjjj] [compare to 

It is important to clearly understand the relation between the stochastic PDE's 
(I30|33l) and the underlying discrete particle systems. The coefficient A in (133]) is defined 
through the current- density relation according to ( [27]) . but the viscosity v and the noise 
strength D in (1291) do not directly appear on the discrete level. To give these coefficients 
a consistent interpretation, we start from the observation [5TL 157] that the invariant 
measure of ( 130]) with spatio-temporal white-noise driving is spatial white noised with 
strength Djlv. This is easy to check for the linearized equation (A = 0) but it remains 
true also for A ^ 0, somewhat analogous to the invariance of the uniform measure for 
the ASEP discussed in Sect. 13.31 As a consequence, the spatial statistics of h(x, t) for 
long times is that of a Wiener process with "diffusion constant" D /Au in spaced , 

limE[(%,£) -h(x',t)) 2 } = —\x-x'\ = A\x - x'\. (34) 

t— s>oo 2l/ 

This relation holds also on the discrete level, provided \x — x'\ is large compared to 
the correlation length of the particle system, and it identifies the ratio A = Dj2v as a 
property of the invariant measure of the latter; for the continuous time ASEP we read 

|| Note that it is possible to have A ^ even if the current itself vanishes at the specific mean density 

under consideration (see [2 3) for an example). 

\ This remains true for certain discretizations of (p?3"|) [751 llllj . 

+ An immediate consequence of (|34l) is that typical configurations of h(x, t) are non-differentiable. This 
is the origin of the ill-posedness of the KPZ equation. 
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off the relation A = p(l — p) from ( 12 5 p . and for the discrete time TASEP A can be 
computed from the transition probabilities ( TTl~|IT2|) [see (|40)) ]. It can be seen from the 
relation ( 1321) [or its discrete analogue ( 12"%|) ] that the height difference correlation function 
is a measure of the fluctuations in the particle number in the interval between x and x' . 
For this reason A has been referred to as a (nonequilibrium) compressibility |54j . 

We note for later reference that the KPZ-equation (1331) can be linearized using the 
Hopf-Cole transformation 



Z(x, t) = exp 



±h{x,t) 



(35) 



which was originally applied to the deterministic Burgers equation [561 [30] and 
rediscovered in the context of (I3"3"j) by Huse, Henley and Fisher [HTJ . Indeed, using (1531) 
we see that Z(x, t) evolves according to a heat equation with a multiplicative stochastic 
force, 



(36) 



The formal solution of ( )36|) is a Wiener path integral describing the weight Z(x,t) of a 
Brownian path (or "directed polymer" [63]) subject to the random space-time potential 



4-2. The universality hypothesis 

The considerations in the preceding subsection suggest that the details of the underlying 
particle system enter the large scale fluctuations properties only through the two 
parameters A and A. These parameters define characteristic scales of height, length 
and time, which can be used to non-dimensionalize any correlation function of interest. 
In the non-dimensional variables the correlation functions are then conjectured to be 
universal, i.e. independent of the specific microscopic model. This is the essence of the 
universality hypothesis. 

As an illustration, consider the probability distribution of the height h(x, t) at a 
given point x, corresponding to the time-integrated current through a fixed bond in the 
exclusion process. Because of translational invariance, this cannot depend on x, and we 
have to find a combination of A, A and t that has the dimension of h. Denoting the 
dimension of a quantity X by [X], we read off from ( 133]) that 

[Al - [X]2 



and from ( 13^1) that 
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The unique combination with the dimension [h] is (A^AIt) 1 / 3 , and hence we expect that 
the rescaled height fluctuation 

k = TOP (37) 
should have a universal distribution. For example, the variance of the height is predicted 
to be of the form [76] 



E[h(x, t) 2 } - {E[h{x, t)}) 2 = c 2 (A 2 |A|t) 2/3 (38) 

with a universal constant c 2 which is independent of the specific model or of model 
parameters such as the update probability it in the dTASEP (see Remark 1 5. 5 1) . Similarly, 
the unique combination of A, A and t that has the dimension [x] of length is 

£(t) = (y4A 2 ) 1/3 t 2/3 , (39) 

which defines the correlation length of fluctuations; note that correlations spread 
superdiffusively, that is, faster than t 1 ^ 2 |137j . This is in contrast to the case A = 0, where 
a straightforward solution of (133]) shows that the correlation length grows diffusively, 
and height fluctuations are Gaussian and of order t 1 ^ [72J. This behavior has been 
explicitly demonstrated for interacting particle systems in which the current is a linear 
(or constant) function of the densitjO [SJ E2 E3]- 

As an illustration of these considerations, and for later reference, we compute the 
scale factor A 2 \X\ for the dTASEP at density p = 1/2. Taking two derivatives of the 
current function (flBT) we find 

2tt 

AdTASEp(l/2) = 



To determine the compressibility A we appeal to the equivalence of the invariant measure 
© to the equilibrium state of the one-dimensional Ising chairj^. Ising spins (Tj are 
canonically related to the occupation variables rji by crj = 1 — 2r]i = ±1. The transition 
probabilities ( fTTTfTST) make up the transfer matrix of the Ising chain, with the density p 
playing the role of the magnetic field (which vanishes when p = 1/2) and the update 
probability 7r controlling the nearest neighbor coupling; since pi/ 2 (0, 1) > pi/ 2 (l, 1) 
the coupling is antiferromagnetic for 7r > 0. Particle number fluctuations translate to 
fluctuations of the magnetization, and hence the compressibility is proportional to the 
magnetic susceptibility of the Ising chain. This can be computed from the free energy 

* It is also possible to construct situations where the leading order nonlinearity in the expansion of the 
current is of cubic or higher order in the density fluctuations, so that A = but the problem remains 
nonlinear [40] . Non-rigorous analysis indicates that such nonlinearities are irrelevant (in the sense that 
the diffusive behavior is preserved) when of quartic order or higher, while in the cubic case fluctuations 
spread weakly superdiffusively as t^Qnt) 1 ' 4 [4T1 IT5] . 

tt The one-dimensional Ising chain is treated in most textbooks on statistical physics, see e.g. [97] . 
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per spin, which is proportional to the logarithm of the largest eigenvalue of the transfer 
matrix, by taking two derivatives with respect to the magnetic field. The final result is 

A dTAS EP(l/2) = \ ^4 = iv 7 ^ (40) 
4 pi/ 2 (0, 1) 4 

and we conclude that 



The early work on KPZ-type processes was mostly concerned with establishing 
the universality of the t 2//3 -scaling of the variance (j38p which, once the role of A and A 
has been recognized, is essentially a consequence of dimensional analysis |72j . Numerical 
evidence of universality in a more refined sense, which encompasses universal amplitudes 
like c 2 in ( 138|) . was presented in [76], where it was also pointed out that different 
universality classes characterized by the same t 2 ^ 3 scaling but different amplitudes 
may arise from different initial and boundary conditions. Specifically, three cases were 
identified: 

I. Growth from a flat surface without fluctuations. In the language of exclusion 
processes, this corresponds to an ordered initial condition of constant density; for 
example, the case p = 1/2 is realized by occupying all odd or all even sites of the 
lattice. 

II. Growth from a flat surface with stationary roughness, in the sense of ( l34l) . This 
corresponds to starting the exclusion process in a configuration generated from the 
invariant measure, e.g. a Bernoulli initial condition of density p for the continuous 
time ASEP. In this case the universal fluctuations of interest are visible only if 
the density is chosen such that the kinematic wave speed c(p) = 0; otherwise they 
will be masked by the fluctuations in the initial condition which drift across the 
observation point. The drift can be eliminated by moving the observation point at 
the kinematic wave speed c(p). 
III. Growth of a cluster from a seed. For the exclusion process this corresponds to a step 
initial condition of the form (120]) with p L > p R . When p L > 1/2 > p R the relation 
(122]) ensures that the density at the origin x = remains at <fi = c _1 (0) = 1/2 at 
all times. As in case II., current fluctuations at other values of p can be studied by 
moving the observation point along a general characteristic x/t = £ with 0(£) = p. 

Early attempts to derive refined universal information, such as amplitudes and 
scaling functions, directly from the KPZ equation met with limited success [76j [2j |52] . 
A full understanding of the universality classes of the one-dimensional KPZ equation 
became available only through the spectacular developments that were triggered a 
decade ago by the paper of Baik, Deift and Johansson [6J. In the next three sections we 
explain the key steps of this development along the lines of the work of Johansson [59], 
and return to the broader issue of KPZ universality in Sect. [SJ 



(^4 |A|)dTASEP 




— 7T. 



(41) 
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5. An exactly solvable model: dTASEP with step initial conditions 

In this section we begin our discussion of Johansson's result [59] on the fluctuations 
of the particle flux for discrete time TASEP (dTASEP) with step initial data. We 
formulate the result in Theorem 15.41 below and we compare it with the predictions of 
KPZ theory described in the previous section. 

Let us first recall the dTASEP model that has been introduced in Sect. 13.11 (iii). 
We denote the infinitely many particles of the system by integers j = 0, 1,2, . . . and 
their respective positions at integer times t — 0, 1, 2, ... by Xj(t) G Z. We assume step 
initial conditions Xj(0) = —j. Jumps to the right Xj(t + 1) = Xj(t) + 1 are attempted 
at every time step t > by all particles j > independently with probability tt, but 
have to be discarded by the exclusion property if at time t the receiving site Xj(t) + 1 
is occupied by another particle of the system, i.e. if = Xj(t) + 1. In this case, 

particle j remains on its site, Xj(t + 1) = Xj(t). 

Definition 5.1 We denote by F n the probability measure on the (total) motion of the 
particle system that is induced by the stochastic process described above. 

Let us first look at an example and compute the probability that the motion 
depicted in Figure [2] occurs. To do this we only need to count for each particle 
j = 0, 1, 2, 3 how many times it had a choice to jump and how often it actually jumped. 





(j choices 


tl jumps 


(j stays 


j = 


10 


5 


5 


J = l 


11 
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i = 2 


10 
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J = 3 
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4 


5 


total 


40 


18 
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By the assumed stochastic independence of all jumps we have 

(Figure [2] occurs) = tt 18 (1 - tt) 22 . 

Next we turn to the flux which is the quantity that we want to analyze. 

Definition 5.2 For r G TL, t G N we denote the total flux through the bond between 
sites r and r + 1 up to time t by 

F r (t) := tJ{j G N: Xj (t) > r} - #{j G N: Xj {0) > r} , 

i. e. the total number of particles that have crossed from site r to r + 1 during the time 
interval [0,t]. 

For example, in the particular situation displayed in Figure [2] we have 
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Figure 2. Sample path for dTASEP where only the motion of the four rightmost 
particles j = 0, 1, 2, 3 is displayed up to some time t < 15. 
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From now on we will only consider the flux F (t) through the bond between sites and 
1 in order to keep the presentation as simple as possible. 

Let us first recall what the discussion on the hydrodynamic limit presented in 
Sect. 13.41 implies for the current at x — 0. We are exactly in the situation of the 
rarefaction wave (see (ii) of Sect. 13.41) with pi, — 1 and Pr = and with J(p) 
given by (USD- Since c(l/2) = J' {1/2) = we learn from (|2l]) and ([22} that 
p(0,t) = 0(0) = c _1 (0) = 1/2 and again by ( fT6l) it follows that the current j (0, t) 
is given by 

j(0,t) = J(p(0,t)) = ^(1 - VT^) =: J.. (42) 
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We therefore expect that F (t) is approximately given by J n t. Indeed, it is a corollary 
of Theorem 15.41 below that F (t)/t converges with probability 1 to J n as t — > oo. 

Now we turn to the more detailed predictions of KPZ theory. As it was explained in 
Sect. I3.5[ the flux F (t) corresponds to the height h(l/2,t). The dimensional analysis of 
Sect. H] suggests that the fluctuations of F (t) — J T t are of order t 1 / 3 (see (1371) ). Moreover, 
taking fl38|) and (TjTJ into account, it is predicted for large values of t that the standard 
deviation of the centered and rescaled variable Z n (t) := t _1 / 3 (F (t) — J^t) is given by 



with C independent of 7T. The result of Johansson forcefully reaffirms the KPZ 
conjecture for dTASEP. Even more is true: In the limit t — > oo, the random variables 
Z^(t) I (7r 1 / 3 (l — 7T) 1 / 6 ) do not only have the same second moments for all < 7r < 1, but 
they converge to exactly the same probability distribution. 

One may think of this result in analogy to the Central Limit Theorem. There one 
considers independent, identically distributed random variables Xj. The quantities for 
which we draw the analogy to the fluxes F (n) are the partial sums S n — Xi + . . . + X n . 
Under some weak assumptions on the distribution of the X^s one has with probability 
1 that S n /n converges to the expectation /z := E(Xl) for n — > oo (law of large numbers) 
and that the rescaled random variables n~ 1 ^ 2 (S n — nfx) tend to a Gaussian distribution 
(Central Limit Theorem). 

In contrast to the Central Limit Theorem the rescaled variables Z v do not converge 
to a Gaussian distribution. Indeed, and quite surprisingly, the limiting distribution is 
given by the Tracy- Widom distribution of random matrix theory. 

Reminder 5.3 (Tracy-Widom distribution) 

The Gaussian Unitary Ensemble GUE is defined as a sequence P^v of Gaussian 
probability measures on N x N Hermitean matrices of the form 



where Zn denotes the norming constant and dM abbreviates the product Lebesgue 
measure on the real and imaginary parts of the entries of M respecting the Hermitean 
symmetry, i.e. dM = UidMu Hi<j dKe(Mij) dlm(Mij) . Denote by Ai(M) the largest 
eigenvalue of M which is a random object. Ai(M) is expected to be at \/2N and 
the fluctuations are of order N' 1 ^ 6 . One can prove that the appropriately rescaled 
largest eigenvalue converges to some distribution that is now called the Tracy-Widom 
distribution, i.e. for s G M we have 



The function TW2 can be expressed in terms of the Hastings-McLeod solution of the 
Painleve II equation |127] or by Fredholm determinants of integral operators with Airy 



V / Var(^(t)) = Ctt 1 / 3 (1 - vr) 1 / 6 , t large 



d¥ N (M) 



—e-'^dM 




as — > 00 . 
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kernel (see Sect. 17.21 for more details). Note that the subindex 2 of TW 2 is related to 
the fact that GUE is a /3-random matrix ensemble with /3 — 2. Roughly speaking a 
/3-ensemble is an ensemble where the joint distribution of eigenvalues is of the form 

i N 
dF N (X u ...,X N ) = -H A(X)f J] w N (Xj)dXj , 
Z N j=1 

and A denotes the Vandermonde determinant (cf. Sect. l6.4l below). See [93] for a general 
reference on Random Matrix Theory. The densities of the Tracy- Widom distributions 
TWi and TW 2 are displayed in Fig. [3]^a) of Sect. EJ 

For our model, dTASEP with step initial data, the theorem of Johansson [59] 
implies: 

Theorem 5.4 Let < tt < 1. Set J : = 1(1 - V^n) and V : = 2- 4 / 3 tt 1 / 3 (1 - tt) 1 / 6 . 
Then, for all s G M we /lave 

^(^I^^)" 1 -™^' < 43 > 

Remark 5.5 ^4s it was noted above Theorem \5.4\ affirms and strengthens the KPZ 
predictions. In particular, the scaling of the flux is precisely that expected from KPZ 
theory. Comparison with (f^7[ ) shows that V = (v4 2 |A|/2) 1 / 3 in the notation of Sect. [^} 

The results of Johansson in [59] are more general than stated above. Mean and 
fluctuations of the particle flux are described not only at the origin and continuous (in 
time) TASEP is also considered by letting n tend to and by rescaling time in an 
appropriate manner. 

We would like to emphasize that Johansson's proof of Theorem 15.41 does not make 
any use of the considerations regarding the hydrodynamic limit and the KPZ conjecture 
as presented above. Instead, the problem treated in Theorem 15.41 should be viewed 
as a very special one within the class of models considered in Sects. |3] and HJ This 
problem has the attractive feature that it is exactly solvable by a series of beautiful and 
non-obvious observations which will be described in the following two sections. 

6. Proof of Theorem 15.41 - part I: Combinatorics 

We begin our discussion of the proof of Theorem I5.4l by relating the flux F (t) to another 
random variable. For j, k e N denote 

T(j, k) := mm{t G N: Xj(t) = k + 1 - j} , 

that is the time by which particle j, that starts at site Xj(0) = —j, has just completed 
its (k + l)-st jump. Observe, that at time := T(k, k) we have 

xo(T fc ) > xi(T k ) > ...> x k (T k ) = 1 > > x k+ i(T k ) > ... . 
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Thus, at time T k exactly the first k + 1 particles 0,1, ... ,k have already jumped from 
site to site 1 and F (Tk) — k + 1. Moreover, for times t < T k we have i*b(^) — k. This 
implies the relation 

K(F (t) < k) — K(T k >t) = l- K(T(k, k)<t). (44) 

In the present section we outline how the explicit formula (1601) of Lemma [6731 below 
for the probability distribution of T(k, k) can be derived. By a series of bijections we 
map our combinatorial model via waiting times and random words to Semi Standard 
Young Tableaux, a classical object of combinatorics and representation theory where 
explicit formulas for counting are available. The asymptotic analysis of formula (|60|) for 
F v (T(k, k) < t) is discussed in Sec. [3 



6.1. From discrete TASEP to waiting times 

We introduce an equivalent description of the dynamics of the particle system by a table 
of waiting times. For j, I 6 N we denote 

Wj t i := number of times particle j decides to stay on site I — j 
after it becomes possible to jump to site I — j + 1. 

For example, in the case of Figure |2] we can determine the following entries of the matrix 
(vjjj) of waiting times. 

\ 



/ o 

2 
1 
2 



3 1 



2 

2 



V 



(45) 



The key observation for computing T(k, k) from the table of waiting times is the 
following recursion for T(j, k). 



T(j, k) = 1 + w j>k + < 



if j = k = 
if j = 0,k > 
if j > 0,k = 
if j, k > 



(46) 





T(j,k-1) 
T(j-l,k) 

max(T(j - 1, k),T(j, k - 1)) 

Indeed, to compute the time it takes the j-th particle to complete its (k + l)-st jump 
one needs to add 1 + Wj tk to the time when this jump became possible. For this jump to 
become possible, particle j has to be on site k — j (happens at T(j, k — 1)) and particle 
j — I must have emptied neighboring site k — j + 1 (happens at time T(j — l,k)). 
It is obvious from (1461) that in order to compute T(k, k) one only needs to know the 
{k + 1) x {k + 1) topleft section of the table of waiting times (ujj,i)o<j,i<k- 
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Relation f !46p allows to prove the following formula for T(j, k) by induction on 

(j + k): 

T(j, k) = j + k + 1 + max I ^ w s ) . (47) 

\sonP / 

Here ILj ^ denotes the set of paths V in the table of waiting times that connect the 
(0, 0)-entry with the (j, fc)-entry and satisfy the additional condition that only steps to 
the right-neighbor and to the neighbor downstairs are permitted. More formally, we 
may write 

U jik = {(s , s j+k ) E (N x NY +k+1 : s = (0, 0), s j+k = {j, k) and 
Si - s^ e {(1, 0), (0, 1)} for all 1 < i < j + k} . 

For V = (so, • • • , Sj + k) G Tlj t k we understand 

j+k 

sonP i=0 

We illustrate formula (T4Tj) with our running example. The corresponding table of waiting 
times displayed in (j4"5l) has nine paths in n 3 3 that maximize the sum of waiting times. 
Two of them are 

P l : (0, 0) (0, 1) (0, 2) (0, 3) -+ (1, 3) (2, 3) (3, 3) (48) 
P 2 : (0, 0) (1, 0) ^ (2, 0) ^ (2, 1) ^ (2, 2) ^ (3, 2) (3, 3) (49) 

and we have 

w s = ^ w s = 8 . 

s on Pi s on V2 

Formula fj47|) then yields for the time when the fourth particle j = 3 has just completed 
its fourth jump T(3, 3) = 3 + 3 + 1 + 8 = 15 which is easily verified from Figure [2j 

Remark 6.1 The probabilistic model we have arrived at, i.e. to search for right- and 
downward paths that maximize the total waiting time, is also known as the last passage 
percolation problem and that is precisely the model studied in the paper 15^ of Johansson. 
Interpreting Wjj as potential energies this can also be considered as the problem of zero- 
temperature directed polymers in a random medium [64\ \7(\ E3 E3, [Tffi . 

We are now ready to compute P 7T (T(k, k) < t) in terms of the table of waiting times. 
Recall from (j4Tj) that T(k, k) is completely determined by the topleft (k + 1) x (k + 1) 
corner of the table. Moreover, and again by (1471) . a table of waiting times corresponds 
to a particle dynamics with T(k, k) < t if and only if the topleft corner belongs to the 
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set W(k, t) which we define to be the set of (k + 1) x {k + 1) matrices (wjj) with entries 
that are non-negative integers and which have the additional property 



max 
vm k . k 



w s) <t-2k-l. (50) 

\s on "P / 

It is straightforward to determine the probability that the topleft corner of the table 
of waiting times agrees with any given element Q of W(k,t). Indeed, we only need to 
count the total number of decisions to either jump (always equals (k + l) 2 ) or to stay 
(equals the sum of all entries of Q which we denote by \Q\i) whenever a jump is not 
prohibited by the exclusion property. In summary we obtain the following formula 

F n (T(k,k) <t) = 7T (fc+1)2 (l-7T) IQ|1 - (51) 

6.2. From waiting times to random words 

We associate with any {k + 1) X (k + 1) matrix Q = (wjj) of waiting times the sequence 
of pairs (j, l)o<j,i<k, listed in lexicographical order, where the value of Wjj determines 
how often the index (j, I) appears in this list. In the case of Q being the top left 4x4 
submatrix in (14"5]) the corresponding sequence of pairs reads 

00001112222233333 
11120030113300223 ( °~' > 
We may consider this list of pairs as a list of 17 two- letter words from the alphabet 
{0, 1, 2, 3} in lexicographical order. This explains the term "random words" often used 
in this context. Observe that any right-downward path V G H33 corresponds to a 
subsequence in this list of 17 two- letter words where both the first and the second row are 
weakly increasing. A little more thought shows that the quantity maxp g n 3 3 (J2 S on v w ^> 
is given by the length of the longest subsequence that is weakly increasing in both rows. 
The sequence (1521) has nine such subsequences of maximal length 8. The subsequences 
corresponding to the paths V\ and V2 of (jUJ) read 



00001223 11222333 
1 1 1 2 3 3 3 3' aU(l 1 1 2 2 3 



Formula (ISTl) translates to 



K(T(k,k)<t)= n^\l-n) len ^ thoi t (53) 

<f>£D(k,t) 

where D(k, t) is the set of finite sequences of lexicographically ordered two-letter words 
from the alphabet {0,1, ... ,k} and for which the length of the longest subsequence of 
that increases weakly in both letters is at most t — 2k — 1. By the Robinson-Schensted- 
Knuth correspondence we may enumerate the set D(k,t) conveniently in terms of Semi 
Standard Young Tableaux. This is the content of the next section. 
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6.3. From random words to Semi Standard Young Tableaux 



The Robinson-Schensted correspondence provides a bijection between permutations and 
Standard Young Tableaux that is well known in combinatorics and in the representation 
theory of the permutation group. We now describe the extension of this algorithm to 
random words which was introduced by Knuth [67] . The basic algorithm that needs to 
be understood first is the row insertion process. Suppose we have a weakly increasing 
sequence of integers, e.g. 1113. We insert an integer r into this row by the 
following set of rules. If r > 3 we simply append r at the end of the row. In the case 
r < 3 we replace the unique number s in the row that is strictly bigger than r such 
that after the replacement the sequence is still weakly increasing. We say that we have 
inserted r by bumping s. For the sequence 1113 insertion of r leads to 



r 


sequence after insertion of r 


bumped number 





000 113 


1 


1 


1111 


3 


2 


1112 


3 


3 


1113 3 


no number bumped 



We now describe the procedure how a lexicographically ordered list of two-letter random 
words is transformed into a pair of tableaux. In a first step, one only considers the 
sequence of the second letters of the random words in order to build the first tableau. 
To this end one inserts (see row insertion process described above) the second letters of 
the words, one after the other, into the first row of the tableau that has been created 
before. In case a number is bumped, the bumped number is inserted into the second 
row of the tableau. In case bumping occurs again in the second row, we insert the 
newly bumped number into the third row of the tableau. This process is repeated until 
the bumping ends. In our example (l52l) this process leads to the following sequence of 
tableaux. 



1 1 



1 1 1 



1112 




1 



1 1 2 



1 

1 1 



00123 o 00023 
11 ^ 1 1 1 



1 3 
1112 



1 1 
1112 3 



1 1 3 
1112 3 



1 1 3 3 
1112 3 



1 3 3 
11113 



Interacting particle systems, KPZ universality, and random matrices 



31 



0000033 0000023 

11111 1 1 1 1 1 3 

2 3 2 3 

0000022 00000223 
1111133-^1111133 

2 3 2 3 

Using this procedure we have obtained a Semi Standard Young Tableau with 17 entries. 



Definition 6.2 By a Semi Standard Young Tableau (SSYT) we understand a tableau 
T of a finite number of integers that are weakly increasing in each row and strictly 
increasing in each column. The shape A = sh(T) ofT is denoted by the sequence of row 
lengths (Xq, X\, . . .) that is required to be a weakly decreasing sequence of non-negative 
integers. Furthermore, we set |A| := t° ^ e the total number of cells in the tableau. 

In order to see that the procedure described above always leads to a SSYT one only 
needs to convince oneself that adding one number to a SSYT will result in a SSYT with 
one more cell. The key observation here is that in the bumping procedure a number can 
only move downwards or left-downwards. 

In our running example the final SSYT T* has shape (8, 7, 2, 0, 0, . . .) and it is no 
coincidence but a theorem that the length of the first row Ao equals the length of the 
longest weakly increasing subsequence that we have seen to be 8. This fact can be 
proven in general by showing inductively that at every step of the procedure the length 
of the longest weakly increasing subsequence (in both letters) up to some word in the 
list is given by the position at which the second letter of the word is inserted in the 
first row of the tableau. In our example (152]) the length of the longest weakly increasing 
subsequence up to the 15-th word (3, 2) is six and correspondingly 2 is inserted at the 
sixth position of the first row. 

Note that the sequence (152]) is not the only one that leads to the final tableau T*. 
For example, the sequence 

00001112222233333 
11123330011300022 

leads to the same T*. However, and this is the central message of the Robinson- 
Schensted-Knuth correspondence, one may encode the sequence of random words ( |52f) 
in a unique way if one records in addition how the tableau grows and if one remembers 
the first letters of the random words that have so far been neglegted. This information 
is all encoded in the second tableau. We now demonstrate how to build this second 
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tableau with our running example. As a first step we record for the above described 
procedure which cell has been added to the SSYT at which step. 

1 2 3 4 7 11 12 17 

5 6 8 9 10 15 16 (54) 
13 14 

Since we also need to remember the first letters of our 17 random words it is natural to 
replace the entries in (l54j) in the following way. We note that the first 4 words in (1521) 
have first letter and we therefore replace 1, 2, 3, 4 each by 0. The next three words 
have first letter 1 and we replace 5, 6, 7 each by 1. Then there are five words starting 
with letter 2, leading us to replace 8, 9, 10, 11, 12 each by 2. The remaining five entries 
13, 14, 15, 16, 17 are each replaced by 3. This leads to the tableau 

1 2 2 3 

1 1 2 2 2 3 3 . (55) 
3 3 

We have arrived at a SSYT U* that clearly has the same shape as T*. In order to 
prove that the second tableau always yields a SSYT there is only one non-obvious 
property to verify, which is the strict increase in each column. To see this recall that 
for words (01,02), {bi,b 2 ) the number b 2 can only bump a 2 if a 2 > b 2 which implies by 
the lexicographical ordering of the list that 61 > a%. Therefore, in the second tableau, 
61 will be located below or left-below a\ and by the weak increase within each row this 
suffices. 

In summary we have mapped the sequence ( 1521) of random words to the pair of 
SSYT's (T*,U*) 

00000223 00001223 
1111133 , 1122233 

2 3 3 3 

of equal shape A. It is an instructive exercise to reconstruct sequence (l52l) from the pair 
of SSYT's. In fact, the proof that the above described procedure maps lexicographically 
ordered two letter words bijectively to pairs of SSYT's of equal shape can be given by 
an explicit description of the inverse map. This bijection has two more features that 
are of interest to us. Firstly, |A| equals the number of words in our list (= 17 in our 
running example). Secondly, the length of a longest weakly increasing subsequence is 
exactly given by the length of the first row Ao (= 8 in our example) as we have already 
observed above (see [67j for details). 

This implies that the set D(k,t) (cf. (1531) ) is bijectively mapped onto the set of all 
pairs (T, U) of SSYT's of equal shape A satisfying 

t - 2k - 1 > A > Ai > . . . 
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with entries from {0,1, ... ,k}. Note that we have A^+i = because entries in each 
column are strictly increasing. We therefore arrive at 

P ir (T(fc, k) <t) — J2 7r (fc+1)2 (l-7r) |A| L(A,fc) 2 , (56) 

t-2fc-l>A >...>A fc >0 

where L(X, k) denotes the number of SSYT's of shape A = (Ao, . . . , A&, 0, . . .) and with 
entries from {0, 1, . . . , k}. We have now derived a representation for P 7r (T(fc, k) < t) 
involving the combinatorial quantity L(X, k) that can be computed explicitly. 

6.4- Schur polynomials and an explicit formula for the distribution ofT(k, k) 

There is a beautiful argument using Schur polynomials s\ that provides a formula for 
L(X,k). This argument is explained in the appendix of [110] and we do not repeat it 
here (see in addition \108\ Cor. 4.6.2] for a derivation of the representation [110[ (A. 5)] 
of s\ by determinants). The result is 

L(\k)= J] K - X i +J -\ (57) 

0<i<j<k J 

Introducing the new variables yi := \— i+k and denoting the Vandermonde determinant 
by A(y) = Uo<i<j<k(yj ~ Vi) we obtain 

k 

K(T(k,k)<t) = C^ k A(y) 2 I](l-7rf\ where (58) 

t-k-l>yo>...>y u >0 i=0 

a,,:=vr( fc+1 ) 2 (l-7r)-^ +1 )/ 2 J] — L- (59) 

0<i<j<k ^ ' 

Observe that each term in the sum is a symmetric function in y that vanishes if two 
components agree. This leads to the final formula in this section for the probability 
distribution of T(k, k) which is closely related to the distribution of the flux -Fo(^) via 

m 



Lemma 6.3 

P„(T{k,k)<t) = j£^ A(y) 2 n(l-T) w • (60) 

0<yi<t-k-l 

This formula should be compared with the formula for the distribution of the largest 
eigenvalue of the Gaussian Unitary Ensemble (cf. Reminder I5.3f) 

N 

Piv(Ai(M) < A) = ^- / A(y) 2 TTe-^ dy (61) 

Z N J(-oo,A] N „'_1 
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with some appropriate norming constant Z^. Observe that this formula has exactly 
the same structure as fl60l) . The role played by the measure e~ x2 dx for GUE is taken 
by the discrete measure YujLoO- ~ 7r )' 7 ^i supported on N for dTASEP. Here Sj denotes 
the ^-distribution concentrated at j. In the next section we discuss how the method 
of orthogonal polynomials can be used to analyze the asymptotics of such types of 
high-dimensional integrals. 

7. Proof of Theorem 15.41 - part II: Asymptotic analysis 

In the previous section we have derived formulae (160]) and (I44p for the probability 
distribution of the flux Fo(t). We now need to analyze this formula asymptotically 
in a regime where t and k both become large and k « Jt + Vt l l z s with s being an 
arbitrary but fixed real number and J, V being defined as in Theorem 15.41 The key to 
the analysis is the observation that the right hand side of ( 1601) is structurally the same 
as the standard formula for the probability distribution of the largest eigenvalue of GUE 
flBTT) and the method of orthogonal polynomials (see Sect. 17. ip can be applied. The role 
played by Hermite polynomials for GUE will be taken by Meixner polynomials in our 
model. In both cases it is the behavior of the orthogonal polynomials of large degree in 
a vicinity of their respective largest zero that matters in the asymptotic analysis. After 
appropriate rescaling this behavior can be described in terms of Airy functions for both 
Hermite and Meixner polynomials (see Sect. 17.21) . On a technical level this explains 
the occurence of the Tracy- Widom distribution TW2 for GUE as well as for dTASEP 
with step initial conditions. We include in Sect. 17.31 a brief discussion of the universal 
behavior of orthogonal polynomials. 

7.1. The method of orthogonal polynomials following an approach of Tracy and Widom 

Almost 50 years ago, Gaudin and Mehta have introduced orthogonal polynomials to 
random matrix theory in order to study local eigenvalue statistics. A very transparent 
version of the now so called method of orthogonal polynomials is due to Tracy and 
Widom [128] and we will briefly outline their approach for our situation. See also the 
recent book of Deift and Gioev [35l Chapt. 4] for a self contained presentation of the 
method of orthogonal polynomials with some remarks on the history of the method. 
For us the method allows to express the right hand side of (1601 in terms of a Fredholm 
determinant of an operator where the (discrete) integral kernel is given in terms of 
Meixner orthogonal polynomials. 

We begin our discussion by introducing the discrete weight 

M*)==f ' • *ez. 

I (1 — 7T) X , if X > 
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For a moment we let (qi)i>o be any sequence of polynomials with qi being of degree / with 
(non-zero) leading coefficient ji. Setting (fi(x) := qi(x) ^ w 7r (x) and using the definition 
of the Vandermonde determinant together with some basic properties of determinants 
we have for y G N fc+1 that 

k 

[det(^))o<,Kd 2 = ( 7o . . . j k ) 2 A(y) 2 J](l - tt) w . 

i=0 

Furthermore, we set I s := [s, oo) and denote by l/ s its characteristic function that takes 
the value 1 on I s and on R \ I s . We may then rewrite ( 1601) in the form 

F n (T(k,k)<t) = - — c ; 2 fc Yl i^^flii-H-M)- 

An identity due to Andreief (see e.g. [351 (3-3)]) adapted to our context reads 

k / 

J2 [det(<^))] [det(^(^))l f[ f{Vi) = (* + 1)! det ( J] ^{x)^{x)f{x) 

One may prove this formula using the Leibniz sum for determinants. This allows to 
write the distribution of T(k, k) as an determinant 

F n (T(k,k)<t)= . Cw ' k . 2 detS, 
(7o • • • 7fc) 

where S denotes the (k + 1) x {k + 1) matrix with entries 

S 3,i = ^2vj( x )Vi( x ) t 1 - l/t-kO*)) > < j,l < k . 

So far the choice of the polynomials qi of degree I was arbitrary. Now we choose (qi)i 
to be the sequence of normalized orthogonal polynomials with respect to the discrete 
measure XLez w tt( x )$x which belongs to the class of Meixner polynomials. We have 

^2<Pj( x )<Pi( x ) = ^2qj( x )qi( x ) w ( x ) = 

x&L xEZ 

for j, I G N. Hence S — I — R(t — k) with 

R ( s )jJ = ^^V^) 1 /*^) = J2vj( x )w( x ) ■ 

x€% x>s 

In summary, we have so far derived 

K(T(k, k)<t)= . C *' k det(/ - R(t - k)). (62) 
(7o • • ■ 7*0 

The prefactor C nt k(lo ■ ■ -1k)~ 2 can be seen to equal 1 by considering equation (1621) for 
fixed k in the limit t — > oo. Thus 

K(T(k,k) < t) =det(I-R(t- k)). (63) 
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The final idea in the argument of Tracy- Widom is to write R(s) - considered as a linear 
map R k+1 ->■ R k+1 - as a product R(s) = A(s)B(s), with 

k 

B{s) : R k+1 £ 2 (Z n J s ), (^)o<i<fc ^ X) WU. 

i=o 

^):f 2 (Zn/ s )^R fc+1 , 

Applying the formula det(J — AB) = det(J — J? A) that holds in great generality 
(see e.g. [351 (3-1)]) we have derived the following Fredholm determinant formula for 
P 7r (T(fc, k) < t). 

Lemma 7.1 P ff (T(Jfe, k) <t) — det(J - S fc (t - fc)), w/jere 

s fc ( s ):£ 2 (zn/ s )^£ 2 (zn/ s ), (X)^?/)/^)) 

and <jfc denotes the reproducing kernel a^{x,y) := £^ =0 V 9 j( x )v 9 j(z/) w ^ respect to the 
Meixner polynomials. 

It may seem somewhat strange to convert (1621) that involves a determinant of some finite 
size matrix I—R into a formula that involves the computation of a Fredholm determinant 
of an operator acting on the infinite dimensional space £ 2 (Z R J s ). However, one has to 
keep in mind that we are interested in an asymptotic result with k — > oo. Hence the 
size of I — R goes to infinity and it is not at all clear how to perform the asymptotic 
analysis of the determinants. In contrast, the operator / — acts on the same space 
£ 2 (Z PI J,,) for all k and it is only the reproducing kernels that dependend on k. As 
will be discussed below the kernels are amenable to asymptotic analysis. In fact, due 
to the Christoffel-Darboux formula for orthogonal polynomials we may express just 
in terms of (fk and (fik+i- For large values of k the behavior of these functions is rather 
well understood. For example, if x is somewhat larger than the largest zero of (pk, then 
is very close to zero. This implies that for values of t — k that are somewhat 
larger than the largest zeros of tpk and (fik+i the operator — k) is negligible and 
thus P 7r (T(fc, k) < t) is very close to 1. If one reduces the value of t — k to lie in a 
vicinity of the largest zero of ifk (which is also close to the largest zero of y^fc+i) then the 
functions <£> fc and (pk+i, appropriately rescaled, are described to leading order by Airy 
functions. In Sect. 17.21 we will use the just mentioned properties of Meixner polynomials 
to complete the proof of Theorem 15.41 

As it was noted in the last paragraph of Sect. 16.41 the formula for the distribution 
of the largest eigenvalue of GUE (16T1) is structurally the same as formula (1601) for the 
distribution of T(k, k) and the arguments described in the present section can be applied 
in an analogous way. The only difference is that we need to use Hermite polynomials 
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instead of Meixner polynomials and that the summation operator is to be replaced 
by an integral operator with a kernel that is given by the reproducing kernel for Hermite 
polynomials up to degree N — 1 (N as in Pjv, cf. Reminder I5.3p . As in the Meixner case, 
the leading order behavior of Hermite polynomials near their largest zero is described 
by Airy functions. On a technical level this is the reason why the fluctuation of the flux 
in dTASEP with step initial conditions follows asymptotically the same distribution 
as the fluctuation of the largest eigenvalue of GUE. It is no coincidence that Meixner 
polynomials and Hermite polynomials of large degree look locally the same when rescaled 
appropriately. In fact, large classes of orthogonal polynomials display the same local 
behavior. We will comment on this universality property of orthogonal polynomials in 
Sect. O 

7.2. Completing the proof of Theorem \5.J\ 

We have argued above that in order to complete the proof of Theorem 15.41 we need to 
determine the asymptotic behavior of the reproducing kernel o~k (see Lemma [7.11) in a 
vicinity of the largest zero of (pk- This requires some detailed analysis that we are not 
going to present here and we refer the reader to (591 Sect. 5]. See also Sect. I7.3l for a few 
general remarks on the asymptotic analysis of orthogonal polynomials. We start with 
some notation. The Airy function can be defined for x G K by 

Ai(x) := — / exp (i[x(t + is) + (t + is) 3 /3])dt (64) 
2vr J R 

with an arbitrary choice of s > 0. The Airy kernel is given by 

x-y 

with the obvious interpretation on the diagonal x = y. The Tracy- Widom distribution 
for (3 = 2 can then be expressed as a Fredholm determinant 

TW 2 (s) := det(I - A)\ L 2 [s>oo) , s£l, 

where A denotes the integral operator associated with the Airy kernel. Note that one 
may derive a differential equation for TW2 leading to another representation |127] 



TW,W=«p(-/"(,-.W^) 



where u denotes the Hastings-McLeod solution of the Painleve II equation u" = 2u 3 + xu 
that is singled out from all solutions of this ordinary differential equation by the 
asymptotic condition u(x) ~ —Ai(x) for x — > 00. Observe that the Airy function solves 

. , . cxp( , -(2/3)x 3 / 2 ) 

the linearized ramleve 11 equation u = xu with asymptotics Ai(x) ~ — 2^7^374 — as 
x — > 00. In an interesting paper [2T] Bornemann demonstrates that it is advantageous 
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for numerical evaluations of 7W 2 to start from the Fredholm determinant formula rather 
than using the Hastings-McLeod function. 

The result in [59, Lemma 3.2] on the reproducing kernel a k for Meixner polynomials 
reads 

ck 1/3 a k (bk + ck 1/3 £, bk + ck 1/3 r]) ->■ A{£, rj) for k ->■ oo, (65) 

where b = 7r _1 (l + yT — n) 2 and c = 7r _1 (l — 7r) 1 / 6 (l + a/1 — 7r) 4//3 . We can now derive 
( H3|) formally. From ([31]) , Lemma 17. 1] and ( 165]) we learn that we should have 

t-fc = 6A; + cA; 1/3 (-s)(l + o(l)) , (66) 

for t — > oo and with k = Jt + Vst 1 ^ 3 . A straight forward calculation shows that this 
can only be achieved if 

(6 + l)J = l and V(b + 1) = cJ 1/3 

leading to the formulae for J and V as presented in Theorem 15.41 Moreover, it is 
apparent that the remainder o(l) in ( 166]) is of order O (t~ 2 / 3 ). Clearly, the argument 
just made does not fully prove Theorem 15.41 since additional estimates are needed to 
deduce the convergence of Fredholm determinants from the convergence of the kernels. 
For details see [59] Sect. 3]. 

Finally we observe that the linear growth of the mean and the t 1 / 3 -scaling of the 
fluctuations of the flux F (t) follow directly from condition (166]) . i.e. from (j4"4"j) . Lemma 
17.11 and ( 16 5 j) without any reference to KPZ theory. 

7. 3. Remarks on the universal behavior of orthogonal polynomials 

We have discussed two examples for fluctuations that can be described by the Tracy- 
Widom distribution: The largest eigenvalue of GUE and the particle flux at the origin 
of dTASEP with step initial conditions. Both systems can be analyzed by the method of 
orthogonal polynomials. The fact that the limiting distribution for the fluctuations agree 
in both cases is then a consequence of the fact that the corresponding sets of orthogonal 
polynomials, Hermite and Meixner, have the same asymptotic behavior in a vicinity of 
their respective largest zeros after appropriate rescaling. This is no coincidence. During 
the past ten years many detailed results on various types of orthogonal polynomials 
have become available that show universal behavior of orthogonal polynomials of large 
degree on a local scale. In this section we give a rough description of this universal 
behavior and outline a few approaches how such results can be proved. 

Let us assume that the support of the measure of orthogonality a is contained in 
R. Then the normalized orthogonal polynomial q n of degree n with positive leading 
coefficient, that is defined through 

J qn{x)q m {x) da(x) = 5 n ^ m , 

k 
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has n simple real roots which we denote by x\ . For many measures of orthogonality a, 
and in particular in the case of the so called varying weights which we do not discuss here 
any further, there exists some natural scaling x — > x such that the counting measures 
n SILi <^-( n ) associated with the rescaled zeros x^ converge for n — > oo to some measure 

i 

H of total mass 1. The support S of /x is compact and is always contained in the support 

2 

of a. For example, in the case of Hermite polynomials da = e~ x dx, the scaling is given 
by x = xjrt}! 2 and the limiting measure of zeros is given by djj, = n~ 1 \ / 2 — x 2 ls(%) dx 
with S =[-V2,V2]. 

One should note that the scaling x — > x as well as the measure of zeros /i do 
depend on the measure of orthogonality a. In order to explain what is universal about 
orthogonal polynomials we restrict ourselves to the common case that S, the support 
of fi, consists of a single interval or a finite union of disjoint intervals. It is convenient 
to describe the behavior of q n by considering the corresponding functions <p n that are 
orthonormal with respect to Lebesgue measure. For example, in the case da = w(x)dx 
we have <p n = q n \fw. In the situation described above the large n behavior of (p n (x) is 
generically the following: 

For x outside S: ip n (x) decays at an exponential rate to zero as n — > oo. 
For x in the interior of S: (p n {x) is oscillating rapidly (in x) and can be described 
to leading order by the cosinus function with slowly varying amplitude and frequency, 
which only depend on fi. 

For x close to a boundary point b of S: To leading order tp n (x) can be expressed in 
terms of special functions. One distinguishes between soft edges (b lies in the interior 
of the support of a) and hard edges {b is also a boundary point of the support of a). 
In the first case the density of n usually vanishes like a square root at b and then the 
leading order of ip n {x) is described by the Airy function. This is in particular the case 
for Hermite and Meixner polynomials. In the case of a hard edge the situation is a bit 
more complicated. In many cases Bessel functions can be used for asymptotic formulas 
for (f n (x) (see e.g. [821 l8T| ). 

We conclude this section by mentioning a few methods how such asymptotics for 
if n can be proved. We will in particular remark on the appearance of the Airy function. 
/. Differential equations of second order 

We again discuss Hermite polynomials da = e~ x dx as a typical example. The 
corresponding functions (p n satisfy the second order differential equations 

^(x) + (2n + l-x 2 K(x) = 0. 

WKB analysis of these differential equations shows that the oscillatory region \x\ < 
y/2n + 1 is connected with the exponentially decaying region |x| > \/2n + 1 by Airy 
functions. This approach can be applied to a number of classical orthogonal polynomials 
that are known to solve linear differential equations of second order with nice coefficients. 
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77. Representation by contour integrals 

Such representations are known for a number of classical orthogonal polynomials 
(e.g. for Meixner polynomials, see e.g. [59j Sect. 5]) and can be analyzed using the 
method of steepest descent. The appearance of the Airy function can be seen from its 
integral representation ( 1641) which generically provides a normal form at critical points 
of higher degeneracy. 

III. Riemann- Hilbert problems 

The characterization of orthogonal polynomials as unique solutions of certain 
matrix Riemann-Hilbert problems (see [33] and references therein) works in principle for 
all types of weights and opens in particular the way to analyze non-classical orthogonal 
polynomials. Here the limiting measure of the zeros fi, which can also be defined as the 
unique minimizer of a variational problem, yields the key to the asymptotic analysis. In 
the neighborhoods of boundary points of the support of /i at which the density vanishes 
like a square root - this is the generic case for a soft edge - Airy functions arise naturally. 
This method for the asymptotic analysis of orthogonal polynomials was first carried out 
in [36]. The method works best in the class of analytic weights, but progress has recently 
been made for weights that have only a finite number of derivatives [91]. Orthogonal 
polynomials with respect to discrete measures have been analyzed by Riemann-Hilbert 
techniques in [SJ. 

IV. Reproducing kernels 

In recent years universality results for orthogonal polynomials, in particular results 
on the reproducing kernel (cf. Lemma [77T1) . have been substantially generalized (see [87] 
and references therein). A very nice view on universality has been introduced in |86] 
where classical results on reproducing kernels for entire functions of exponential type 
are being used (see [M] for the Airy kernel). 

8. KPZ-universality revisited 

In view of the universality conjecture formulated in Sect. HJ one expects that the results 
derived in the preceding sections for a very special case - the dTASEP with step initial 
conditions - should carry over, in a quantitative sense, to a much broader class of models. 
The first explicit demonstration of this idea was presented by Prahofer and Spohn in 
a series of papers [99, I100[ 1101] , where an alternative and independent route linking 
Ulam's problem to growth models was established 

The starting point is the one-dimensional polynuclear growth model (PNG), an 
interacting particle system on the real line, in which particles (antiparticles) move 
deterministically at unit speed to the right (left), annihilate upon colliding, and are 
created in pairs according to a two-dimensional Poisson process in space and time [78J. 

f Another link between the two classes of problems was found by Majumdar and Nechaev [89l 188] . 
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Figure 3. (a) The densities of the three universal distribution functions TW2, TW\ 
and Fo (from left to right). Discrete points show simulation results for the PNG model. 
Reprinted with permission from jlOOj . Copyright 2000 by the American Physical 
Society, (b) Phase diagram for the distribution of current fluctuations in the TASEP 
with Bernoulli step initial conditions. Here G denotes the Gaussian distribution and 
G 2 (GOE 2 ) is the distribution of the maximum of two independent Gaussian {TW\) 
random variables. Initial particle densities are ph — P- to the left and pu = p+ to the 
right of the origin. Reprinted from [101) with kind permission of Springer Science and 
Business Media. 

Via the random set of particle creation events the model can be mapped onto the 
problem of the longest increasing subsequence of a random permutation, which in turn 
provides a link to the Tracy- Widom distribution [6]. For the case of a droplet growing 
from a seed (case III of Sect. H]), Prahofer and Spohn show that the resulting fluctuation 
distribution is identical (under the rescaling prescribed by KPZ theory) to that obtained 
by Johansson for the dTASEP. 

Moreover, by imposing suitable boundary conditions [9] and symmetry relations 
[TO] on the set of Poisson points, the cases of flat and rough initial conditions (case 
I and II of Sect. |4|) can be handled as well [100J. For the flat initial condition (case 
I) the fluctuations are governed by the GOE distribution TWi, while for the rough 
initial condition (case II) a new distribution Fo emerges which so far does not have 
an interpretation in terms of random matrix theory [9]. The three distributions are 
depicted in Fig. Eta). 

Since the fundamental works of Johansson, Prahofer and Spohn the field has 
developed rapidly, and it is impossible to do justice to the new results in the framework 
of these lecture notes. In the following subsections we therefore restrict ourselves to 
briefly outlining the most important directions of research, providing the interested 
reader with a few key references along which recent advances can be traced. 
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8.1. The Prahofer-Spohn conjecture and the ASEP 

Based on the universality hypothesis, Prahofer and Spohn translated the results obtained 
for the PNG model into a conjecture for the fluctuations of the particle current through 
the origin for the TASEP with a general step initial condition (120]) . where particles are 
placed to the left (right) of the origin according to a Bernoulli measure with density 
Pl (pr) [101j . The fluctuation phase diagram in the plane of the boundary densities 
is shown in Fig. ]3](b). The overall features of the diagram can be understood from 
hydrodynamics. First, the Johansson result obtained at p^ = 1, pr = is seen to extend 
throughout the region p L > 1/2 > p R . As explained in Sect. HI this reflects the fact that 
the density profile near the origin is independent of the boundary densities in this case. 
For pi < 1/2 and pn > 1/2 the application of the hydrodynamic formulae (1211) [for 
Pl > Pr] and ( |T9i) [for pi < pr] show that the density at the origin becomes pi and pr, 
respectively. In these cases the intrinsic current fluctuations are masked by the initial 
fluctuations drifting across the origin, leading to simple Gaussian statistics (regions 
marked G in the diagram). The line Pr + Pl = 1, Pl < Pr, is special, because there the 
shock speed f|T9|) vanishes and the density at the origin shifts randomly between p L and 
Pr. As a consequence, the current is distributed as the maximum of two independent 
Gaussian random variables (denoted by G 2 in the figure). Similarly, along the lines 
Pl = 1/2, pr < 1/2, and pr = 1/2, pi > 1/2, the distribution is that of the maximum 
of two independent variables drawn from TW\. Finally, at the point p^ = pr = 1/2 we 
have case II behavior governed by the distribution Fq. 

A proof of the Prahofer-Spohn conjecture for the TASEP was recently presented by 
Ben Arous and Corwin [4] (see also [961 El SS] for earlier partial results). Moreover, in a 
remarkable series of papers Tracy and Widom have been able to generalize these results 
to the (partially) asymmetric exclusion process [131] 11301 11321 1135] 11331 1134] 1136] . The 
generalization is highly nontrivial, because the ASEP for general q is not a determinental 
process [133], and it requires a novel set of techniques based on the Bethe ansatz 

[TrgJ E05J [53]. 

8.2. Spatio-temporal scaling 

We have seen in Sect. HJthat the essence of the KPZ conjecture is the universality of 
height fluctuations when viewed on the appropriate scales defined by the height rescaling 
f]37]) and the correlation length f]39]) . In other words, once the average growth shape has 
been subtracted, one expects that the rescaled fluctuations 

h t {y) = {A 2 \\\t)-y*h{y{A\H 2 ) l '\t) 



converge for t — > oo to a universal stochastic process A(y), whose single-point 
distribution is one of the random matrix distributions discussed above. The process 
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A(y) was first explicitly characterized by Prahofer and Spohn for the PNG-model in the 
droplet geometry (case III of Sect.H]), who named it the Airy process |102j . Subsequently 
Sasamoto identified the analogous process for the case of flat initial conditions (case I) 
[109J. In line with the nomenclature used to designate the corresponding single-point 
distributions (TW 2 for case III and TW\ for case I), the two processes are now called 
Airy 2 and Airy! processes, respectively [2JJ [25J [231 EQ] • Whereas the Airy 2 process has 
a natural interpretation in the random matrix context as the motion of the largest 
eigenvalue in GUE matrix diffusion (Dyson's Brownian motion), the corresponding 
relation does not hold for the Airy! process [22]. The corresponding process for initial 
conditions with stationary roughness (case II) was studied in [7] (see also [58]). 

In recent work Ferrari and collaborators have extended the analysis to include 
correlations between height fluctuations at different times t and t! > t. As for the 
scaling of the height fluctuations themselves, the characteristics of the hydrodynamic 
equation play a special role for the decay of correlations. Whereas the decay along 
generic space-time directions is governed by the correlation length ( 1391) . which is of 
order t 2//3 , along the characteristics the decorrelation time at time t is set by t itself, 
which implies a much slower decay |l9j [32j [31]. This is in accordance with KPZ 
phenomenology, which predicts that such correlations should decay as {t/t') x with a 
universal autocorrelation exponent A [691 EIJ- In contrast to the scaling exponents 
of single-point height fluctuations introduced in Sect. HI the autocorrelation exponent 
depends explicitly on the growth geometry: For a flat initial condition (case I) A = 1, 
whereas for a curved cluster (case III) A = 1/3 |121j . 

8.3. KPZ-scaling at large 

The results described so far in this section were based on a small set of exactly solvable 
models, the (T)ASEP's and the PNG model. On the other hand, KPZ universality is 
expected to hold for a much broader class of interacting particle systems and growth 
models which is limited only by the requirement of local, stochastic transition rules and 
a nonlinear dependence of the particle current (or growth rate) on the particle density 
(or surface slope) [see Sect.H]. It is therefore gratifying that the class of models for which 
KPZ universality has been rigorously established - mostly in the sense of finding the 
exact scaling exponents governing the order of fluctuations - has been greatly expanded 
in recent years. This has required the development of new techniques that are purely 
probabilistic in nature and do not rely on the specific analytic structure of the exactly 
solvable models. The first result of this type was obtained in [28] for the Hammersley 
process, an exlusion-type process in continuous space which is closely related to the 
Ulam problem. Similar methods were subsequently applied to a variety of interacting 



Interacting particle systems, KPZ universality, and random matrices 44 

particle systems |1U3[I1U4| [TS | [TB]. including a class of zero range processed with general 
jump rates that have a non-decreasing, concave dependence on the number of particles 
|12j . This constitutes a major step on the way to proving KPZ universality in the 
broadest sense. 

8.4- The universality class of the KPZ equation 

Ironically, although tremendous advances in the analysis of different representatives of 
the KPZ universality class were achieved over the past decade, the one-dimensional KPZ 
equation ( 133]) itself remained rather poorly understood. We noted already that the KPZ 
equation is mathematically ill-posed because of the highly singular white noise term, 
and some regularization is needed to make it amenable to rigorous analysis. This can 
be done by spatial discretization or by constructing the equation through a scaling limit 
from an asymmetric exclusion process with weak asymmetry [T7] 1112] . Both approaches 
have recently been used to prove the correct order of fluctuations in the KPZ equation 
[111] [13] as well as refined universality in the sense of the Tracy- Widom distribution 
|1131 11141 [3]. Thus it has finally been established, as it were, that the KPZ equation 
belongs to its own universality class. 

In another line of work an independent, non-rigorous approach to establishing 
Tracy- Widom universality has been developed which is based on applying the replica 
method to the path weight Z(x,t) in the Hopf-Cole transformed equation (136|) . In this 
approach one computes moments E(Z n ) with respect to the stochastic force ((x,t), 
which results in a problem of n bosonic, quantum-mechanical particles interacting 
through an attractive ^-function potential [621 [55] . While the ground-state energy and 
wave function of this quantum system been known for a long time, the recent works 
[271 H2] have succeeded in summing over the full spectrum of excited states of the 
many-body Hamiltonian. 

8.5. Experiments 

Despite the wide applicability of the KPZ theory to a large class of stochastic growth 
models, experimental signatures of KPZ scaling in the real physical world have proven 
to be surprisingly elusive. Early efforts focused on the investigation of growth-induced 
(two-dimensional) surface roughness of crystals and thin solid films [70J. However, 
detailed consideration of the physical processes governing such growth experiments has 
revealed that they typically operate in regimes where KPZ asymptotics is practically out 
of reach [721 OS]- To date, the most thoroughly studied experimental system that displays 
KPZ scaling is the slow, flameless combustion of paper. By imaging and analyzing the 



| Zero range processes were introduced by Spitzer |122j and have been extensively studied in the physics 
literature, see [44] . 
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Figure 4. (Color online) Experimental demonstration of KPZ universality in a thin 
film of turbulent liquid crystal, (a) Outlines of a growing turbulent droplet at time 
intervals of 5 s. (b) Measured distribution of shape fluctuations in comparison to 
the GUE prediction. Reprinted with permission from [126] . Copyright 2010 by the 
American Physical Society. 

one- dimensional smoldering front, the t 1//3 -scaling of fluctuations |90j as well as the non- 
Gaussian shape of the corresponding probability distribution [95J were demonstrated. 

Very recently, a new experimental system involving different phases of driven, 
turbulent liquid crystal films became available, in which the refined universality 
predictions of the theory can be tested with unprecedented precision |126j . In Fig. H] we 
show a series of snapshots of the growing turbulent droplet along with the experimentally 
determined distribution of shape fluctuations. Under the rescaling prescribed by KPZ 
theory, the predicted universal GUE distribution is seen to emerge from the data without 
any adjustable parameters. 

9. Integrability and Universality 

The central model discussed in our paper is the totally asymmetric simple exclusion 
process (TASEP). This model has two features that motivate our choice. Firstly, it 
belongs to the class of stochastic interacting particle systems that are useful in the 
study of transport phenomena in nonequilibrium systems. We are particularly interested 
in the description of fluctuations around the mean behavior which is governed by 
a deterministic evolution equation. In general it is difficult to obtain such detailed 
information, mainly because the interactions between the particles destroy stochastic 
independence and the classical central limit theorem cannot be applied. However, and 
this is the second feature, TASEP is an exactly solvable model and the fluctuations 



Interacting particle systems, KPZ universality, and random matrices 



46 



can be analyzed by a series of beautiful and non-obvious observations. We have seen 
by explicit calculation that the fluctuations of the flux are described by a formula that 
is structurally the same as the formula for the fluctuations of the largest eigenvalue of 
matrices from the Gaussian Unitary Ensemble (GUE). This provides one explanation 
for the much celebrated link to Random Matrix Theory. There are also links to the 
theory of integrable systems (see e.g. [33], I139[ 1140] ) that we have not explained at all 
and that we mention here only in order to give an indication why TASEP is sometimes 
called an integrable model, even though there is no differential equation to integrate in 
sight. 

A common feature of integrable systems is that their delicate mathematical 
structure immediately breaks down if the model is changed ever so slightly. One may 
think that this limits the interest in the corresponding results. However, the recurring 
experience with integrable systems has been, that even though the method of proof is 
not applicable if the model is changed, the results may persist. In our context a nice 
example of this principle of universality is the recent work of Tracy and Widom where 
it is shown by very different methods that a number of results for TASEP also hold for 
ASEP. Universality results are also available for large classes of matrix ensembles (see 
[351113] and references therein). Despite these results, the question of universality is a 
subtle one as can be seen for example from the correlations of the fluctuations (cf. [22] 
and Sect. 18.2ft . It remains a great challenge to understand the realm of validity of the 
various laws that have been established for the specific integrable models. 
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